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A metastable lattice gas with nearest-neighbor interactions and continuous-time dynamics is stud- 
ied using a generalized Becker-Doring approach in the multidimensional space of cluster configu- 
rations. The pre-exponential of the metastable state lifetime (inverse of nucleation rate) is found 
to exhibit distinct peaks at integer values of the inverse supersaturation. Peaks are unobservable 
(infinitely narrow) in the strict limit T — » 0, but become detectable and eventually dominate at 
higher temperatures. 
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In a general case, the nucleation-controUcd lifetime of 
a metastable state can be written as 

; T ^ Aexp{W-,/T) (1) 

with VF* being the minimal work (free energy change) 
to form a critical nucleus and the temperature T mea- 
sured in units of Boltzmann constant. The exponential 
term was anticipated already in the earliest estimations 
[01 of the nucleation rate / ~ l/r, which further found 
an enormous amount of applications in systems ranging 

' from vapors [|| to glass-forming [|j or quantum |^ liq- 
uids. The structure of the pre-exponential, A, however, is 
not known in a general case. The difficulty of its deriva- 
tion, whether in an analytical, numerical or experimental 
study, stems from the dominant contribution of the ex- 

[ ponential in eq.(^, with minor uncertainties in (say, 
due to inaccuracies in the measured interfacial tension) 

I implying large, orders of magnitude, discrepancies in the 
values of A. 

With this, much attention is devoted to models which 
can exhibit nucleation and which are close to exact solv- 
ability. Here one can obtain an accurate expression for 
W*, subsequently focusing on the pre-factor issue. 

One of the best known example is the two-dimensional 

[ nearest neighbor Ising model, where metastability is 
achieved by orienting initially all spins one way (down) 
while non-zero magnetic field h prescribes an oppo- 
site (upward) orientation. Allowing spin flips of non- 
conserved Q or conserved |^ type adds the required dy- 
namics to the problem. In a closely related lattice gas 
model the role of h is played by supersaturation. 

The pre-exponential in such systems attracts much 
attention both for the high- ||, ^, |l^ and low- 
temperature [jll], ^ regions. For h the nucleus 
is macroscopic and its shape, as well as the value of W^, 

' can be obtained from the Wulff droplet construction . 
Monte Carlo simulations are possible for < 10 — 15 T 
which, for small h, restricts such studies to the afore- 
mentioned high-temperature region; transfer-matrix ap- 
proaches also are available here For larger fields a 
straightforward Wulff construction may be inadequate, 
Ip^ , but for T — > analytical treatment becomes possi- 
ble due to dominant contribution of low-energy configura- 
tions. Technique of absorbing Markov chains [0 [l2[ |l^ 



also can be used for simulations in the low-temperature 
region. 

Neves and Schonmann[T7t evaluated W^, the zero- 
temperature limit of W*, obtaining the exponential part 
of the metastable lifetime. Their result is insensitive to 
specifics of the dynamics. Novotny further showed 
that for discrete-time dynamics and a relatively large 
field, the pre-exponential remains finite in the limit T — > 
0, approaching a piece- wise constant function of h, point- 
ing towards a discontinuity at an integer value of the in- 
verse field. Similar features will be observed at weaker 
fields as well Integer values of inverse field, however, 
were excluded from the aforementioned rigorous mathe- 
matical treatments, leaving open questions with regard 
to this intriguing effect, especially in the physically more 
realistic case of T > 0. 

The present Letter aims to evaluate the pre- 
exponential at higher temperatures and in a finite domain 
of fields, spanning several integer values of 1/h. This 
will clarify the nature of the discontinuities and, together 
with the available , will provide predictive expressions 
for t{T, h) at T > 0. We will show that in contrast to 
intuitive expectation of discontinuities spreading out in 
a standard, tanh-\ike fashion, they are replaced by sharp 
peaks which persist, with finite heights and self-similar 
shapes, up to T = 0. 

A lattice gas model with continuous time dynamics 
will be considered. Specifically, the probability of cre- 
ation of a particle on an empty site in an infinitesimal 
time interval dt is taken as pdt, regardless of the sur- 
rounding; without restrictions, the time scale I3~^ can 
be taken as 1. Alternatively, the annihilation probabili- 
ties are proportional to dtex-p{—AE/T), with AE being 
the energy change due to broken bonds and field ("su- 
persaturation" ) which increases the energy by 2h when a 
particle is removed. This model, with various generaliza- 
tions, is popular, e.g. in Monte Carlo simulations of the 
dynamic interface in crystallization problems - see, e.g. 
|18[ and references therein. It is expected that qualita- 
tively the model also remains similar to the discrete-time 
Glauber type dynamics of Refs.jl^, |l^ (and the general- 
ized Becker-Doring approach employed below bears cer- 
tain parallels with the technique of absorbing Markovs 
chains ||ll|, ^), although the dynamics-sensitive pre- 
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exponentials will not be identical even at T = 0. 

For a long time of the order of t the rare particles 
will form isolated clusters of various sizes and shapes 
{classes), which will be distinguished by a running in- 
dex, i. An empty site corresponds to i — 0. Cluster 
shapes will be considered identical (and thus belonging 
to the same class) if they can made such by rotation or 
reflection. The key characteristics of each class are the 
numbers of particles, s{i), the number of bonds b{i) and 
the statistical weight Wi < 8. One can define the (quasi) 
equilibrium distribution 

f^l = u;,z2s(0~6(«)js(^) (2) 

with z = e~'^^'^ and S = z~^^ describing the temperature 
and field dependencies, respectively, (p is the bond en- 
ergy, subsequently taken as 1 for simplicity of notations. 
In the s, b space the function /f' has a saddle point (for 
non-special fields - a single one ) and the correspond- 
ing value of s determines the critical cluster number, s*. 
In a general case computer assistance is required in order 
to characterize all classes. Consistency of such predic- 
tions can be checked, e.g., against standard tables 
for smaller s. 

Once equilibrium properties are specified, one can in- 
troduce kinetic fluxes as a multidimensional version of 
the classical approach since in a low-temperature 
Ising system growth or decay of a cluster predominantly 
proceeds via random gain or loss of a single particle [gO| . 
If Pikdt is the probability to transform a cluster from 
class i to class k > i by adding a particle [with (3ik = 
if s{k) ^ s{i) + I], the corresponding flux is given by 

hk = ftfc/r(Wi -Vk) , i < k (3) 

with Vi = fi/ fl'^ and wq = 1. The Master Equation for 
the kinetic distributions fi takes the form 

f = E^-- E (4) 

fc=0 k=i+l 

which automatically satisfies detailed balance. 

For closing conditions, absorbing states are placed at 
all classes k with s{k) — Smax + 1- Equivalently, all those 
absorbing states can be combined in a single absorbing 

class ^niax f ■ 

Due to an exponentially long lifetime, one can neglect 
the depletion of empty sites (for which, otherwise, an 
integral conservation law should be employed instead 
of Wo = !)■ With this, eqs. (|), (g) can be solved in the 
steady-state approximation; transient effects |Q also can 
be neglected here. 

Introducing hk = Pikff^ + Pkifl'^ (0 < i, fc < fcmax + 1) 
and 

femax+l 

Mik = bik - Sik ^ bil , 1 < Z, fc < fcinax (5) 



one can show that the steady-state distributions 
(vi , V2, . . . ) correspond to the first column of the matrix 
—M~^ (since only class 1, with single-particle clusters, 
is connected to empty sites). The total flux I coincides 
with Iqi where branching of paths does not yet occur. 
This gives 

/=(M-i)n + l (6) 

For a single nucleation path, which leads to a tri- 
diagonal structure of the matrix M, one recovers the 
classical result by Farkas / ^ = 691^ + ^12 + • ■ • • 

Oth- 
erwise, the actual evaluation of / via eq. (^) is limited 
by one's ability to obtain all classes and transition rates 
f3ik for a sufficiently large Smax, and the ability to in- 
verse analytically a large matrix M . At present, we were 
able to proceed up to Smax = 9 (1818 classes represent- 
ing a total of 13702 cluster configurations) which allows 
us to consider fields h > 1/6 with the critical number 
s* < 7. A full exact expression for r can be surveyed by 
a human eye only for more modest values of Smax, which 
implies a relatively small critical cluster (larger fields). 
For example, for s,nax = 4 kinetics is determined by 9 
distinct classes with a total of 28 shapes (see, e.g. Fig. 
2 in Ref.[|o|). Transition rates are easy to obtain (say, 
there arc two ways a 3-particle "minus" shaped cluster 
can turn into 4-particle "T" shaped one, four ways it can 
turn into an "L" shaped one, etc.). The result which fol- 
lows from eq.(^), is expressed as a rational function of z 
and 5 

Ti^PiS,z)/Q{S,z) (7) 

with the subscript indicating the value of Smax, and poly- 
nomials P and Q given by 

P{d,z) = 384 + 210000(5^z^ -I- 165(48 + 185z) + (8) 
4:6^z{1576 + 2655z) + 86^z^{308l + 3230z) + 
25005^^^^(21 + 250z + 362^) + 
2505^2*^(695 + 2650z + 10362^) + 
5^^3(65740 + 57797Z + 15360^^) + 
55^z^(28574 + 28155z + 19680^^) -f 
55^^(43375 + 70546z -t- 50000^^) 

g(5,z) = 85S*(384-|-8O(24-|-855)0-h 

2505^(25 -I- 645)z3 + 1255^(259 + 6256)z'^ + 
205(615 + 17535)^2 + 37505'*(3 + 75)z^) 

Eq.(^ is expected to be accurate in strong fields, 
h > 1/2, with rather relaxed restrictions on tempera- 
ture since all cluster configurations at s < 4 are taken 
into account (although, for higher T eventual destruc- 
tion of the steady-state due to neglected cluster interac- 
tions should be kept in mind |23|| ). More consistently, 
this result should be treated asymptotically for z ^ 
and 5 — » 00 with certain combinations of powers of z and 
5 remaining finite, depending on the interval of field. 

In order to isolate the pre-exponential, eq.(0) should 
be multiplied by exp(— W,/T) = z^* . In principle, an 
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"observable" is r itself, rather than A or W^, taken sep- 
arately. To avoid ambiguity, the value of will be 
associated with its zero-temperature limit, 0, with 
all temperature-dependent corrections being in the pre- 
exponential; for h > 1 the barrier will be taken as zero. 
The function W^(h) has a piece- wise linear structure, 
and exp(— W*/r) is reduced to a product of integer pow- 
ers of z and 5: 1 for /i > 1, z^S for 1/2 < /i < 1, z'^S^ 
for 1/4 < h < 1/2, etc. The resulting A{h) is shown by 
a dashed line in Fig. Q where numerical results, given 
as filled circles, were obtained for a much larger Smax 
and can be treated as "exact" in the present context. 
The case T = would correspond to a piece-wise con- 
stant structure of A, similar to Ref. jl^ but with differ- 
ent constants and an additional "excluded singularity" 
at h= 1/2: A = 1 for /i > 1, A = 1/4 for 1 > > 1/4 
{h ^ 1/2), .4 = 1/16 for 1/4 > h > 1/6, etc. (These 
numbers can be deduced from the lowest energy path 
-see below-, serving as a checkpoint for more more elabo- 
rate expressions). This limit, however, becomes apparent 
only at a very low temperature, z = 10~^ . 

For a larger cut-off, simplifications of analytics can be 
achieved due to the dominant contribution of low-energy 
configurations. Among all classes k of clusters with the 
same s{k), one can select only those which have a suf- 
ficiently large number of bonds: bs{k) > bmax,s{k) ~ 
where hmax,s is the number of bonds in the most com- 
pact cluster for a given s. An integer parameter r indi- 
cates how close a cluster should be to the most compact 
configuration in order to be included in the kinetics. For 
sufficiently large r (r = 4 for s„iax = 9) all configura- 
tions are recovered. Alternatively, r — Q corresponds to 
the lowest energy path description, which is the closest to 
the kinetic part of the conventional one-dimensional ran- 
dom walk approach to nucleation Q , although with mi- 
croscopic rather than phenomenological coefficients. In 
addition, branching of paths is added starting from s ~ 7 . 
Already in the r — Q approximation peaks at integer l/2h 
will appear in the pre-exponential, although one needs to 
include r > 1 for correct evaluation of their heights. 

For the case r = 1 and z <C 1 the pre-exponential 
A(h) can be described analytically in restricted domains 
of fields, the most interesting being those near the peaks 
(general expressions for A{h) are also available, but are 
useless due to their size). 

Introducing a finite combination 

y = 5z''^ (9) 

with n = 1,2,... determining a corresponding peak, one 
can perform analytical expansions of l/I in fractional 
powers of z. Symbolic computations with Mathematica 
were used here. 

For n ^ 2 one obtains 

1 T^jy) 1 T2{y) 1 T^{y) 

^ z^r2 8y'^T{y) z^ Wy^T^y) z^^ Q72y^T^{y) 

(10) 



with the coefficients in this Smax = 9 approximation given 
by 

T(y) = 8 + 63y2 (U) 
Ti{y) = 4 + 2 V + 79?;^ -M26/ 
T2{y) = -208 - 1432y2 + 3461y'' 49855/ 

89649?/^ + 87318?;i" 
T3,{y) = 89152 - 297840^^ - 13174644?;'* - 62801445?/^ + 

146767614?/^ -f 1284356493yi° + 

9576800102/12 -I- 5566046222/1"* 

The approximation works accurately in the vicinity of 
h = 1/4, describing the rather complex near- and off- 
peak behavior - see Fig. ^j. The coefficient of z^^/^ in 
eq.(|l^) multiphed, respectively by y"^ at ft, > 1/4 or by j/*" 
at /i < 1/4, determines the scaling structure of the peak 
in the limit T — > 0, if the difference /i— 1/4 scales together 
with temperature. The structure of the neighboring peak 
at h = 1/2 {n — 1) follows the 4-particle approximation, 
eq.(|^), with 5 ~ y/z and z — )■ 0. 

An important question is the sensitivity of the results 
to variations in r and Smax- A smaller r — Q can re- 
produce the coefficient of z^^/^ in eq.(p^) in the lim- 
its j/ ^ oo or ?/ ^ (i.e. on both sides of the peak 
for T 0), but will not give the proper peak height 
for y — 1, or scaling for finite y, or the correct higher- 
order terms. Cases with larger r > 2 presently could 
be studied only numerically and are shown by symbols 
in Fig. ^. At small h scatter appears in data, indicat- 
ing the limits of numerical accuracy for very small z. 
There is no detectable difference with the analytical ap- 
proximations in the regions of their validity for fields up 
to /i > 1/4. An similar expansion in z for r = 1 and 
Smax = 8 also was performed, leading to a rather differ- 
ent structure of the y-dependent polynomials. The first 
coefficients in the z-expansion are nevertheless numeri- 
cally close for stronger fields h> 1/4 in the vicinity of the 
peak. So are the heights of the peak given, respectively 
by 0.4167 - 2.92zi/2 -t- . . . and 0.4190 - 2.832^/2 + ... 
in the 8-and 9-particle approximations. On the other 
hand, unlike the 9-particle case, Smax = 8 does not yield 
a proper T — > limit for weaker fields /i < 1/4 since the 
boundary here is too close to the critical size s* — 7. 

In summary, for a moderate field (supersaturation) the 
metastable state lifetime of a supersaturated lattice gas 
has been evaluated for T <^ Tc. The main result is 
the pre-exponential which, for the first time, was evalu- 
ated analytically beyond the zero-temperature limit, and 
which exhibits distinct peaks as a function of field. One 
can anticipate that similar peaks (which appear due to 
competition of several "critical sizes") also will be ob- 
served in systems other than nearest-neighbor Ising mod- 
els with non-conserved dynamics, whenever the nucle- 
ation barrier has a well-defined zero-temperature limit 
and the critical nucleus contains a reasonably small num- 
ber of particles. 
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FIG. 1: The pre-exponential of eq. (|l|). Lines - the 4-particle 
approximation, eq.(^), at a higher temperature, z — 10~^ 
(short dashed) and at a low temperature, z = 10~^ (long 
dashed). Points - numerical results for Smax ~ 9 and r = 2. 
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FIG. 2: "Fine structure" of the two peaks near h — 1/2 and 
h — 1/4 for z = 10"''. Solid line- the 9-particle approximation 
for r = 1, eq.(|l^). Dashed line- the full 4-particle approxi- 
mation, eq.(^). Symbols - numerical data for Smax = 9 and 
r = 2 (small circles) and r = 4 (large circles) , respectively. 
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